Multiple paternity, fertilization success, and male quality: Mating system variation in the eelgrass, Zostera marina

Abstract Genetic diversity can modulate a population's response to a changing environment and plays a critical role in its ecological function. While multiple processes act to maintain genetic diversity, sexual reproduction remains the primary driving force. Eelgrass (Zostera marina) is an important habitat‐forming species found in temperate coastal ecosystems across the globe. Recent increases in sea surface temperatures have resulted in shifts to a mixed‐annual life‐history strategy (i.e., displaying characteristics of both annual and perennial meadows) at its southern edge‐of‐range. Given that mating systems are intimately linked to standing levels of genetic variation, understanding the scope of sexual reproduction can illuminate the processes that shape genetic diversity. To characterize edge‐of‐range eelgrass mating systems, developing seeds on flowering Z. marina shoots were genotyped from three meadows in Topsail, North Carolina. In all meadows, levels of multiple mating were high, with shoots pollinated by an average of eight sires (range: 3–16). The number of fertilized seeds (i.e., reproductive success) varied significantly across sires (range: 1–25) and was positively correlated with both individual heterozygosity and self‐fertilization. Outcrossing rates were high (approx. 70%) and varied across spathes. No clones were detected, and kinship among sampled flowering shoots was low, supporting observed patterns of reproductive output. Given the role that genetic diversity plays in enhancing resistance to and resilience from ecological disturbance, disentangling the links between life history, sexual reproduction, and genetic variation will aid in informing the management and conservation of this key foundation species.


| INTRODUC TI ON
Genetic diversity is a fundamental characteristic of populations, determining their response to a changing environment (Fisher, 1930) and playing a critical role in ecological function (Hughes et al., 2008;Whitham et al., 2006).Especially in systems where an organism serves as a foundation species, the functional traits displayed by different genotypes can have profound effects on community persistence and productivity (Reusch, 2006;Reynolds et al., 2012;Wimp et al., 2004) in ways comparable to species diversity in other systems (Crutsinger et al., 2006;Wimp et al., 2004).While multiple processes act to maintain genetic diversity, sexual reproduction remains the driving force, particularly over short time scales, and is highly dependent upon life-history strategy and mating system (Williams, 1975).
Seagrasses (i.e., marine flowering plants) act as foundation organisms for highly productive nearshore habitats and provide essential ecosystem services such as habitat for an array of epifaunal species and fisheries, coastal protection from storm surge and wave action, sediment stabilization via root and rhizome growth, and the mitigation of excess nutrient loads and eutrophication (Duarte et al., 2013;Gillanders, 2006).Seagrass meadows are, however, experiencing a trajectory of decline, and a global net loss of 5602 km 2 (19.1% of surveyed meadow area) has occurred since 1880 (Dunic et al., 2021).Rates of seagrass loss are comparable to those reported for mangroves, coral reefs, and tropical rainforests, which places seagrass meadows among the most threatened ecosystems on Earth (Waycott et al., 2009).On average, only 37% of restoration efforts are successful (Van Katwijk et al., 2016), highlighting the importance of conserving existing meadows before collapse occurs.
The seagrass Zostera marina (eelgrass) is commonly found in the temperate regions of the world (Green & Short, 2003) and experiences a range of abiotic and biotic stressors (Kendrick et al., 2019;Lefcheck et al., 2018;Orth et al., 2017).Such selective pressures can lead to changes in life-history strategy and mating system dynamics (Cabaço & Santos, 2012).As such, Z. marina's allocation to sexual reproduction varies across its geographic range (Phillips et al., 1983).Specifically, extreme environmental conditions result in annual populations characterized by yearly shoot mortality, development of only reproductive shoots, high seed production, and biennial establishment of seedlings (Phillips & Backman, 1983); within-range, perennial populations are less reliant on flowering and sexual reproduction, additionally undergoing clonal extension of the rhizome to form ramets (Silberhorn et al., 1983).However, even in the same region, flowering can vary widely within and among populations for reasons still largely unknown (Von Staats et al., 2021).
Perennial populations maintained by asexual, clonal growth were once assumed to be the primary life history strategy of Z. marina (Den Hartog, 1970).Now, sexual recruitment from seed is becoming more widely acknowledged as an important driver in both the maintenance and expansion of perennial (Johnson et al., 2020) and annual populations (Xu et al., 2018).
At its southern edge-of-range, Z. marina populations experience annual loss of biomass due to heat stress and are shifting primarily to a mixed-annual life-history strategy characterized by annual, complete loss of biomass, development of both vegetative and reproductive shoots, reestablishment from seeds alone, and seedlings that flower in their first year of growth (Allcock et al., 2022;Bartenfelder et al., 2022;Jarvis et al., 2012).Because they rely heavily on sexual reproduction, southern edge-of-range meadows may represent valuable genetic diversity hotspots among Z. marina populations (Diekmann & Serrão, 2012).Importantly, in Z. marina, genotypic diversity has been associated with resistance and resilience to natural disturbances (Hughes & Stachowicz, 2004, 2009, 2011;Reusch et al., 2005), and allelic richness and relatedness have both been linked to increased biomass and productivity (Hughes et al., 2016;Stachowicz et al., 2013).Therefore, the mating system can be especially influential in genetic diversity hotspots by impacting standing levels of variation within these populations.
Measures of the genetic mating system include polyandry (i.e., the degree to which females multiply mate), paternity skew (i.e., the extent to which mating can be monopolized by a fraction of the available males), and outcrossing rates (i.e., the proportion of crossfertilized versus self-fertilized progeny).Currently, we have a limited understanding of basic characteristics of the sexual mating system in Z. marina.Of the relatively few eelgrass mating system studies, the majority have investigated population-level patterns in clonal structure and outcrossing rates.Results reveal population-specific patterns in clone structure (Billingham et al., 2007;Furman et al., 2015;Peterson et al., 2013;Reusch, 2001), related to environmental impacts on life-history strategies (Harwell & Rhode, 2007).Additionally, there are population-specific patterns in outcrossing rates in which meadows range from almost entirely outcrossing (Furman et al., 2015;Reusch, 2000b;Ruckelshaus, 1995) to nearly entirely selfing (i.e., selffertilizing) in instances of low genotypic diversity (Reusch, 2001).
Furthermore, some populations contain spatially organized kin groups (Billingham et al., 2007;Furman et al., 2015;Hays et al., 2021), suggesting that kin interactions have the potential to shape seagrass behavior, life-history strategy, and resource allocation.
Fine-scale analyses on individual-level patterns in mating system have been conducted in perennial meadows and report multiple paternities within inflorescences, confirming that marine angiosperms are polyandrous (Follett et al., 2019;Hays et al., 2021;Reusch, 2000b).Several studies have also shown that withinmeadow hydrodynamics play a significant role in shaping mating systems of Z. marina of the Northwestern Atlantic, in which outcrossing rates increase in deep water and near the top of the meadow canopy (Follett et al., 2019, Hays et al., 2021).Importantly, the presence of multiple and varying paternities within a maternal brood (in this case, an inflorescent spathe) can reflect potential differences in male siring success and levels of inbreeding.For example, Z. marina genet size in the German Baltic Coast was positively correlated with heterozygosity (Hämmerli & Reusch, 2003a).Similarly, self-fertilization has been shown to decrease seed set in Z. marina (Hämmerli & Reusch, 2003b), suggesting that the impacts of mating system can extend to other plant traits, such as seed size and germination rates, which results in important ecological ramifications for meadow health (Delefosse et al., 2016).
Because edge-of-range populations experience annual stress and display a shifting reliance on sexual reproduction (Jarvis et al., 2012), analysis of the mating system as a driver of genetic diversity is needed.Here, we investigate within-plant patterns in polyandry, siring success and outcrossing in three mixed-annual eelgrass meadows by genotyping all seeds within a given reproductive shoot.This approach overcomes the limitations of previous work by greatly increasing the sample size per reproductive shoot, enabling an assessment of both structural and temporal patterns of mating system variation.By illuminating the links between life history, sexual reproduction, and genetic variation, this study represents an in-depth characterization of within-plant genetic mating system dynamics of Z. marina at its southern limit along the western North Atlantic Ocean.More broadly, as eelgrass systems globally have flowering densities and seed outputs that are similar in magnitude despite their differing life-history strategies (Combs et al., 2021), our results may thus be representative across a wide range of meadow types and locations.
Spathes contain a spadix with rows of both male and female flowers clustered in ratios of 2:1.Flowering occurs in two stages: (1) the flowering and exposure of ovaries and (2) the opening of anthers and release of pollen, each separated by approximately 5 days (De Cock, 1980).Past studies describe sequential flowering of spathes (i.e., temporal dichogamy) within Z. marina reproductive shoots beginning with basal structures (De Cock, 1980, 1981).Specifically, one spathe per rhipidium flowers at a time, during which the next spathe of the same rhipidium is not yet entirely developed.The second spathe starts flowering a few days following the first spathe opening its anthers.Temporal dichogamy seemingly does not exist across rhipidia (De Cock, 1980, J. Jarvis, unpubl. data).In the present study, spathe position is therefore used as a proxy for reproductive timing, and rhipidium position is used as a proxy for height within the water column, consistent with previous literature (e.g., Follett et al., 2019;Furman et al., 2015).

| Study sites and sample collection
Three seagrass meadows located in Topsail Sound, North Carolina, were sampled to characterize Z. marina mating system variation.
The meadows were in a shallow coastal lagoon (34.22 N, 77.37 W) protected from the Atlantic Ocean by Topsail Island (Figure 1).This North Carolina Wilmington's Center for Marine Science.Shoots were then cleaned and blotted dry, and the morphological position of each seed was recorded.Each seed was given a position within a spathe, assigned to a given spathe and to a given rhipidium.Labels were assigned in ascending numeric order from increasing proximity to the rhizome (e.g., basal positions were given a value of 1) (Figure 2).
Seeds were removed from spathes, cleaned, blotted dry, and tested for viability using the "squeeze test" by gently compressing individual seeds with a pair of tweezers (Marion & Orth, 2010).
Those with a seed coat that compressed were considered nonviable.Viable seeds were removed from the seed coat and genotyped.
Historically, the "squeeze test" is performed on seeds that have been released from the flowering shoot; however, if this had been done, the position of the seed within and the identity of the parent plant would not have been known.Moreover, seeds were sampled in May at the peak of the eelgrass reproductive season in North Carolina where most sampled seeds were mature (i.e., at stages 4 and 5; Combs et al., 2021).Indeed, when the sites were re-visited one-week post-sampling, all seeds had been released from their spadices (Jarvis, pers. obs.) indicating that seeds collected as part of this study were mature at the time of collection.Moreover, any immature seeds released from the spadix within 1 week were effectively nonviable; they would not have successfully germinated.

| Seed genotyping
DNA was extracted from viable seed samples using a PowerPlant® Pro DNA Isolation Kit.Ten microsatellite loci previously characterized for Z. marina (Oetjen et al., 2010;Oetjen & Reusch, 2007;Reusch et al., 1999; Table A1) were amplified in two multiplex Polymerase Chain Reactions (PCR).Individual primer working stocks contained 1 μL of 10 μM fluorescently labeled forward primer and 10 μL each of 50 μM unlabeled forward and reverse primers diluted in 80 μL of ddH 2 O. Primers were then combined into two primer mixes -each containing five different primers (Table A1).PCR conditions for all multiplex conditions were as follows: 95.0°C for 15 min; 2 cycles of 94.0°C for 15 s, 60.0°C for 30 s, 72.0°C for 45 s; 2 cycles of 94.0°C for 15 s, 59.0°C for 30 s, 72.0°C for 45 s; 2 cycles of 94.0°C for 15 s, 58.0°C for 30 s, 72.0°C for 45 s; 2 cycles of 84.0°C for 15 s, 57.0°C for 30 s, 72.0°C for 45 s; 28 cycles of 94.0°C for 15 s, 56.0°C for 30 s, 72.0°C for 45 s; and a final 2 min extension at 72.0°C.Following PCR, two reactions were prepared: one containing 0.5 μL of each PCR product from each of the multiplex mixes.PCR products were added to 9 μL of highly deionized formamide (HiDi) and 0.4 μL of GeneScan-600 (LIZ) size standard (Applied Biosystems, Foster City, CA, USA) for capillary sequencing on an ABI Prism 3130XL Genetic Analyzer.Fragments were scored using Applied Biosystems Microsatellite Analysis Software (ThermoFisher Scientific Inc.).

| Paternity analyses
Known maternal half-siblings were used for sibship reconstruction and paternity assignment with the maximum-likelihood approach in COLONY v2.0.7.0 (Jones & Wang, 2010).COLONY parameters included a polygamous mating system for both sexes, inbreeding, and a monecious, diploid species.A long run with medium-likelihood precision and a genotyping error rate of 1% was performed.Maternal and paternal genotypes were reconstructed using an allele probability threshold of 0.925 for allele calls at each locus.Seeds were categorized as selfed if the putative father had the same genotype as the putative mother.Outcrossing F I G U R E 2 Labeling system of Zostera marina flowering shoots.Rhipidia (shown in blue, teal, pink, and green colors) form branching structures on each flowering shoot.Spathes (i.e., inflorescences) of each rhipidium (shown only for Rhipidium 2 in varying shades of green) constitute a branch of a rhipidium and contain a spadix with flowers and eventually seeds.Seeds 1-n rates were calculated as the proportion of outcrossed offspring per meadow, shoot, rhipidium, and spathe.The effective number of sires and paternity skew per spathe, rhipidium, and shoot were calculated after Neff et al. (2008) in which effective sires = 1/Σ(rs i / seeds) 2 where rs i = the number of offspring assigned to sire i, and the summation is over all sires contributing to a maternal spathe, rhipidium, or reproductive shoot.Skew was then expressed as 1 -(effective number of sires/actual number of sires).As such, a value of 0 implies no skew in which case all sires contribute equally to a seed set, and a value approaching 1 implies maximal skew in which case nearly all offspring are assigned to a single father.The paternity skew of each sire was calculated as the proportion of genotyped seeds per shoot sired by a particular sire.

| Statistical analyses
Statistical analyses were conducted in RStudio with R v4.2.1 (Posit Team, 2023; R Core Team, 2022), and figures were generated with "ggplot2" and "lattice" (Sarkar, 2008;Wickham et al., 2016).Data were tested for outliers, collinearity, even sample size, and normal distribution (Zuur et al., 2007).To assess patterns in seed viability, a generalized linear mixed model (GLMM) was fit to test the fixed effects of meadow, rhipidium position, and spathe position on the number of viable seeds per spathe (Poisson distribution, offset by the number of seeds per spathe).A random effect of maternal identity was added to account for differences among mothers.To assess patterns in mating system variation, GLMMs with the appropriate distribution were fit to test the fixed effects of meadow, rhipidium position, and spathe position on the response variables of number of outcrossed seeds (Poisson distribution), paternity skew (log + 1 transformed, Gaussian distribution), and number of sires (Poisson distribution) per spathe using the package "lme4" (Bates et al., 2015).
A random effect of seeds per spathe was added to control for expected variance due to the fair raffle process in sperm competition (Parker, 1990;Zuur et al., 2007), and random effect of maternal identity was added to account for differences among mothers.Model residuals were visually inspected for normality and homogeneity using the package "DHARMa" (Hartig & Lohse, 2022).Global models were used to perform ANOVAs (α = .05)and post hoc pairwise comparisons (α = .05)using the package "car" (Fox & Weisberg, 2011).
To explore whether the number of sires was influenced by (a) the number of genotyped seeds, (b) the total number of seeds, and (c) the proportion of selfed offspring, GLMMs (Poisson distribution) were also fit on both a "per spathe" and "per shoot" basis, with a random effect of maternal identity.To assess the impact of paternal genotype on reproductive success, GLMMs (Poisson distribution) were fit to test paternal heterozygosity, whether a sire selfed or outcrossed, and their interaction on (a) the number of seeds sired per male and (b) paternity skew per male, with a random effect of the number of successfully reconstructed loci per sire.Model residuals were visually inspected for normality and homogeneity using the package "DHARMa" (Hartig & Lohse, 2022).
Overall, 93% of offspring loci were successfully amplified and scored (7849 of 8440 loci).All loci were polymorphic, ranging from two alleles (ZM10, meadow 3) to nine alleles (ZM3, meadows 1 and 2) (Table A3).Given the small quantity of DNA extracted per seed (~5 ng/μL in a 10 μL volume), re-runs were rarely possible.The analyses presented below were performed on this complete dataset of 844 seeds.COLONY reconstructed 251 of 260 maternal loci (97%) with an average probability of 0.997 and 1322 of 2080 paternal loci (64%) with an average probability of 0.994.Additional analyses on a reduced data set of 688 seeds that were successfully genotyped at ≥9 loci were also performed and gave qualitatively similar results (see Tables A1-A5 for analyses on the reduced dataset).
The same patterns were observed on a "per shoot" basis (Table 4).

| DISCUSS ION
Zostera marina meadows growing at the southern edge-of-range in the western North Atlantic experience a high degree of thermal stress and annual summer mortality (Bartenfelder et al., 2022;Thayer et al., 1984).As a result, recent shifts from a primarily perennial to a mixed-annual life-history strategy have been observed (Jarvis et al., 2012), along with resultant increases in flowering densities (Combs et al., 2021;Jarvis et al., 2012;Thayer et al., 1984) and genetic diversity (Allcock et al., 2022).Our characterization of the mating system in three regional meadows revealed high levels of multiple paternity and outcrossing rates across individual reproductive shoots.Mating systems also varied temporally, with some spathe positions showing significantly higher levels of selfing.
Moreover, fertilization was not evenly distributed among individuals as reproductive success varied significantly among males -those with higher heterozygosity and those that self-fertilized produced greater numbers of seeds and displayed higher reproductive skew (i.e., monopolized a greater fraction of the available seeds).Our findings reveal an edge-of-range eelgrass mating system consistent with shifts in life-history strategy toward greater sexual recruitment and highlight the impact of male traits on reproductive output.
Multiple mating occurred in all sampled shoots and was not significantly different among the three Z.marina meadows.

TA B L E 3
Global generalized linear mixed model results for the effects of meadow, rhipidium position, and spathe position on mating system characteristics.
angiosperms (Pannell & Labouche, 2013), particularly in dense populations (e.g., Friedman & Barrett, 2011).For example, up to nine sires within a single fruit have been reported in the terrestrial angiosperms scarlet gilia (Ipomopsis aggregata) (Campbell, 1998), white campion (Silene latifolia) (Teixeira & Bernasconi, 2007), and Allegheny monkeyflower (Mimulus rigens) (Mitchell et al., 2005).A fair raffle process in sperm competition (Parker, 1990) predicts a positive correlation between the number of seeds and number of sires.That is, larger reproductive shoots are more likely to contain offspring from each of the males a female mates with, simply because of the chance probability that each male's pollen contributes to at least one of the offspring.Indeed, this is reflected in the positive correlation detected between both the number of seeds (total and genotyped) and number of sires.Polyandry was also influenced by the proportion of selfed seeds (Table 4).In instances of high degrees of self-fertilization, few sires monopolize available ovules, resulting in lower levels of polyandry.Therefore, the observed variance in the number of sires appears to be a function of the fair raffle process and the degree to which the parent plant selfed.In general, polyandrous mating systems may be advantageous as a mechanism for producing genetically diverse progeny (Karron & Marshall, 1990) and selecting for male phenotypic traits that provide the greatest offspring fitness (Lankinen et al., 2009;Schlichting et al., 1990) which, in turn, may lead to population growth (Ashman et al., 2004).
Paternity skew varied markedly within and among flowering shoots and was not correlated with meadow, rhipidium, or spathe position.Some sires fertilized up to 2400% more offspring than others (Figure 3), and such variance in reproductive success is often attributed to variance in either offspring mortality or male quality (Hedgecock, 1994;Hedrick, 2005).Differential patterns of offspring mortality could contribute to paternity skew among males, with some males producing fewer viable seeds than others.As genotyping was only performed on viable seeds that were close to maturity at the time of sampling, the paternal influence on nonviable seeds remains unknown.While there are reviews of hydrophilous (i.e., water-mediated) pollination in the seagrasses (e.g., Ackerman, 2000;Furman et al., 2015;Ruckelshaus, 1996), there have been no attempts to link pollination to male genotype.Male fertilization success in angiosperms is presumed to be impacted by a variety of factors, including pollen production (Field et al., 2012), pollen clumping (Martin et al., 2009), and pollen competition (Mulcahy & Mulcahy, 1975).
Indeed, intraspecific variation in pollen production and viability has been observed in terrestrial angiosperms and linked to male genotype (i.e., some males produce more or better pollen than others), resulting in some males contributing significantly more offspring than others (Danti et al., 2022;Devasirvatham et al., 2012).Additionally, pollen clumping can increase fertilization success in wind-pollinated species -clumped pollen fertilizes nearby conspecifics, while single pollen fertilizes plants at a greater distance (Martin et al., 2009).
Seagrasses are hydrophilous pollinators and the formation of pollen clouds, similar to pollen clumps, has been suggested in Z. marina populations from Shinnecock Bay, New York (Furman et al., 2015) and in other seagrass species (McConchie & Knox, 1989).Because sires that self-fertilized (regardless of heterozygosity) had consistently higher fertilization success than outcrossing males, proximity between anther and stigma also appears to significantly increase males' fertilization success.While the aforementioned pollen features have not been directly linked to male genotype in Z. marina, the positive correlation between male heterozygosity and fertilization success suggests that male quality may play a role in the observed patterns of reproductive variance.Although heterozygosity is not uniformly found to be a reliable proxy for reproductive success (Botero-Delgadillo et al., 2020;Wetzel et al., 2012), in seagrasses, it has been linked to increased genet size (Hämmerli & Reusch, 2003a).Larger genets, in turn, possess more flowering shoots and thus more pollen, though estimates of clone size and clonal range have not been made in our system.
When Zostera in the western Baltic Sea self-fertilized, the theoretical fitness of selfed offspring decreased (Reusch, 2001); and in the Ria Formosa, relatedness of parent plants increased rates of seed abortion (Billingham et al., 2007;Zipperle et al., 2011).Alternatively, selfed mating produced more viable seeds in the Chesapeake Bay (Rhode & Duffy, 2004), an area where eelgrass experiences summer heat stress (Hensel et al., 2023) and is recovering as water clarity improves (Lefcheck et al., 2018).
Outcrossing rates did not differ across rhipidium positions, which are often used as proxies for canopy height.In contrast, previous studies in Z. marina found that the topmost rhipidium of each shoot displayed high levels of outcrossing, as pollen dispersal was aided by increased water movement at shallower depths within the water column (Follett et al., 2019).Canopy-mediated hydrodynamics had no observed effect on Z. marina mating systems in Topsail, NC, likely due to differences in canopy height and meadow submergence between these meadows and those in northern Massachusetts (Follett et al., 2019).Meadows in northern Massachusetts had a canopy height of 85-135 cm and were submerged in a tidal height of 3 MLLW (Follett et al., 2019); meadows in NC had a canopy height of 5-35 cm and were submerged in a tidal height of 2 MLLW (Jarvis et al., 2022).The shorter canopy height likely reduced the effect of hydrodynamics on mating system variation, and the shallower water depth lessened the effect of above-canopy water velocity.Hays et al. (2021) also found that meadow depth itself was a significant driver of mating system dynamics.Shallow depth (1 m MLLW) reduced water movement and pollen dispersal, resulting in decreased outcrossing and paternal richness compared to deeper sites (3 and 5 m MLLW).The physical structure of seagrass meadows serves to increase water velocity above the canopy and decrease water velocity within the canopy, and the effect of meadows on water velocity is exaggerated in deeper meadows (Fonseca et al., 1983;Thayer et al., 1984).Therefore, pollen disperses farther when more deeply submerged (Fonseca et al., 1983, Thayer et al., 1984).
Outcrossing rates did differ among spathes, suggesting that mating system characteristics in Z. marina can change over time, even within a single plant.Given that selfing appears to increase in younger spathes, within-shoot phenology may provide a potential explanation for this pattern.When the first spathes open on a given flowering shoot, pollen must come from another ramet.
However, as subsequent spathes mature and extend their stigmas, previous spathes are now extending their anthers, allowing for within-shoot pollen transfer (De Cock, 1980, 1981).Additionally, temporal changes in pollen availability can occur through quantity (e.g., not enough pollen to go around) or quality (e.g., decreased viability, inbreeding depression) limitation (Aizen & Harder, 2007).
Pollen quantity limitation has been reported in the Zosteraceae in the Netherlands (Van Tussenbroek et al., 2016) and in the Baltic Sea (Reusch, 2003).The greater proportion of nonviable seeds on younger spathes observed in this study might thus be the result of decreased pollen availability as the season progresses, though viability cannot unequivocally be distinguished from immaturity here as seeds were removed directly from the spathe rather than collected post release.
Seagrass meadows in NC experience frequent natural disturbance from tropical storms (Paerl et al., 2019;Zhang et al., 2022), thermal stress (Bartenfelder et al., 2022;Jarvis et al., 2012), and sedimentation (Mills & Fonseca, 2003), which are known to promote increases in sexual reproduction (Cabaço & Santos, 2012).The high levels of multiple paternity and outcrossing rates are consistent with meadows experiencing high levels of flowering.Moreover, the observed effects of male heterozygosity and temporal shifts in both selfing and viability highlight the contribution of fine-scale processes to population-level patterns of mating system variation, genetic diversity, and, ultimately, ecological function (Wimp et al., 2004).Recent literature has emphasized the need to conserve existing meadows before collapse occurs (Van Katwijk et al., 2016) and to prioritize genetically diverse meadows for their enhanced ecosystem services (Reynolds et al., 2012).Z. marina edge-of-range meadows of the western North Atlantic may represent one such priority, and our findings reveal processes underpinning the genetic composition of this species, which is essential for informed conservation and management efforts.writing -review and editing (equal).

ACK N OWLED G EM ENTS
The authors would like to thank the University of North Carolina Wilmington for financial support.Special thanks to Randy Taylor, UNCW's Coastal Plant Ecology Lab, and the 2022 cohort of UNCW's MarineQuest Oceans 17 for field and laboratory assistance.With respect to the land and people who preceded us, the coauthors of this article acknowledge that this research was conducted on colonized land, traditionally and ancestrally belonging to Native Peoples, all of whom have stewarded this land throughout the generations.We also acknowledge this coastal region profited from the enslavement and forced labor of African people and their descendants.We recognize that knowing, acknowledging, and honoring the history of the land and the people is only the first step and that the continuation of addressing policies and practices that perpetuate oppression is needed.
A PPEN D I X TA B L E A 1 Microsatellite forward and reverse primer sequences.Note: Significant predictors (p < .05)are indicated in bold.Dataset 1 used reconstructed paternal genotypes from the COLONY run with all seeds (n = 844); here only sires with an average genotype probability of ≥.925 were used.Datasets 2 and 3 used reconstructed paternal genotypes from the COLONY run with seeds for which ≥9 loci successfully amplified (n = 688).Dataset 2 used only paternal loci that were reconstructed with a probability of ≥.925.Dataset 3 used sires with an average genotype probability of ≥.925.

Microsatellite
lagoon contains Z. marina in both monospecific and mixed assemblages (i.e., co-occurring with Halodule wrightii and Ruppia maritima) (NCDEQ, 2021).Sites were classified as shallow subtidal (depth <2 m MLLW) and were on a narrow shelf between the Intracoastal Waterway and the adjacent shoreline.Twelve flowering shoots were haphazardly collected roughly 5 m apart from each site (n = 36 shoots total) on May 4, 2021, and transported on ice to the University of F I G U R E 1 Location of Zostera marina meadows sampled within the Topsail Sound Intracoastal Waterway, North Carolina.

F
I G U R E 4 Mean proportion of viable seeds per spathe as a function of spathe position.A, AB, B, and C refer to group bins based on post hoc within-group pairwise comparisons (α = .05).Individual spathe values are plotted using the jitter function in R. F I G U R E 5 Mean proportion of outcrossed seeds per spathe as a function of spathe position.A, AB, and B refer to group bins based on post hoc within-group pairwise comparisons (α = .05).Individual spathe values are plotted using the jitter function in R.

TA B L E 5
Generalized linear mixed models results for the effects of paternal heterozygosity on the number of seeds fertilized and paternity skew.F I G U R E 6 Relationship between paternal heterozygosity and (a) the number of seeds fertilized and (b) paternity skew.Triangles and dashed line correspond to selfing sires; circles and solid line correspond to outcrossing sires.

Table 3 )
, decreasing from 0.86 ± 0.03 at spathe position 1 to 0.52 ± 0.07 at spathe position 4 Variables include the number of genotyped morphological structures (n), the proportion of outcrossed seeds per genotyped seeds x ± SE, the number of sires (x ± SE), the number of structures fertilized per sire (x ± SE), and paternity skew (x ± SE).
Note: Variables include n (sample size), k (mean kinship coefficient), H O (observed heterozygosity), H E (expected heterozygosity), and N A (number of alleles).
Generalized linear mixed models results for the effects of the total number of seeds, the number of genotyped seeds, and the proportion of selfed seeds per spathe and per shoot on the number of sires detected.Significant predictors (p < .05)are indicated in bold.
Sample sizes for each reproductive shoot used in the reported dataset (n = 844).Global generalized linear mixed model results for the effects of paternal heterozygosity and selfing on the number of seeds fertilized and paternity skew using various datasets.